%
% Script que realiza la simulacion para el inciso a.
% Parametros de entrada:
%
%   tfin:        duracion de la simulacion
%
%   ks_func:     nombre de una funcion que retorne un vector de valores
%                de k para repetir la simulacion
%
%   postfijo:    un string que puede ser '_derivativo', 
%                '_integrador' o '_proporcional' 
%
%   ivp: ivp que se desea simular
%
% No retorna nada, genera los archivos de datos para ser graficados
% por a_graph.m
%

function a(tfin, ks_func, postfijo, ivp)

	t0 = 0;
	if (ivp == 'ivpi')
		ci = [0 0 0];
	else
		ci = [0 0 ];
	end

	% Parametros del modelo
	b = 0.02;
	myks = feval(ks_func);
	i = 0.004;

	prefijo_theta = strcat('../datos/theta', postfijo);
	prefijo_tiempo = strcat('../datos/tiempo', postfijo);
	prefijo_theta_r = strcat('../datos/theta_r', postfijo);

	% Parametros para la integracion numerica.
	h = 0.01;
	% Simular para cada k y plotear los resultados
	for j = 1:length(myks)
		j
		params = [b myks(j) i];
	
		% Simular el modelo para los distintos k
		[ts xs theta_r_dat] = rk4(ivp, params, t0, tfin, ci, h);
		thetas = xs(:, 1);
		% Guardar los resultados de la simulacion
		aux = strcat(prefijo_theta, num2str(myks(j)));
		save(aux, 'thetas');
		aux = strcat(prefijo_tiempo, num2str(myks(j))); 
		save(aux, 'ts');
		aux = strcat(prefijo_theta_r, num2str(myks(j))); 
		save(aux, 'theta_r_dat');
	end
	% Simular el blanco
	clear ts;
	clear thetas;
	ts = 0:h:tfin;
	len = length(ts);
	for j = 1:len
		t = ts(j);
		thetas(j) = feval('theta_r', t);
	end

	save('../datos/theta_r', 'thetas');
	save('../datos/tiempo_r', 'ts');
	
end
